%% Simulate Simple Shock Responses
% by Jaromir Benes
%
% Simulate a simple shock both as deviations from control and in full
% levels, and report the simulation results.

%% Clear Workspace
%
% Clear workspace, close all graphics figures, clear command window, and
% check the IRIS version.

clear;
close all;
clc;
irisrequired 20140315;

%% Load Solved Model Object
%
% Load the solved model object built in `read_model`. Run `read_model` at
% least once before running this m-file.

load read_model.mat m;

%% Define Dates
%
% Define the start and end dates as plain numbered periods here.

startDate = 1;
endDate = 40;

% ...
%
% Alternatively, use the IRIS functions `yy`, `hh`, `qq`, `bb`, or
% `mm` to create and use proper dates (with yearly, half-yearly, quarterly,
% bi-monthly, or monthly frequency, respectively).
%
%    startdate = qq(2010,1);
%    enddate = startdate + 39;

%% Simulate Consumption Demand Shock
%
% Simulate the shock as deviations from control (e.g. from the steady
% state or balanced-growth path). To this end, set the option
% `'deviation='` to true. Both the input and output database are then
% interpreted as deviations from control:
%
% * the deviations for linearised variables are defined as $x_t -
% x_t$: hence, 0 means the variable is on its steady state.
% * the deviations for log-linearised variables are defined as $x_t / \Bar
% x_t$: hence, 1 means the variable is on its steady state, or 1.05 means
% it is 5 % above it.
%
% The function `zerodb` automatically detects the maximum lag in the model,
% and creates the input database accordingly so that it includes all
% necessary initial conditions.

d = zerodb(m,startDate:endDate);
d.Ey(startDate) = log(1.01);
s = simulate(m,d,1:40,'deviation=',true);
s = dbextend(d,s);
s %#ok<NOPTS>

s1 = simulate(m,d,1:40,'deviation=',true);

%% Report Simulation Results
%
% Use two quick-report functions, `dbplot` and `qplot`, to report
% simulation results. The latter relies on a q-file, i.e. `shock_report.q`,
% which defines what and how to display in the figure; see help on
% `qreportlang` for more details. The two functions produce identical
% graphs except for titles (use the option `'titles='` in `dbplot` to
% modify the titles the function prints).

plotrng = startDate-1 : startDate+15;
series = {'100*(dP^4-1)','100*(R^4-1)','100*(Y-1)', ...
   '100*(N-1)','100*(W/P-1)','100*(Pk-1)'};
dbplot(s,plotrng,series,'tight',true);
ftitle('Consumption Demand Shock');

qplot('simulate_simple_shock.q',s,plotrng,'tight',true);
grfun.ftitle('Consumption Demand Shock -- Deviations from Control');

%% Simulate Shock in Full Levels
%
% Instead of deviations from control, simulate now the same shocks in full
% levels. To that end, create an input dabase with the steady state
% (balanced-growth path) using `sstatedb`, and keep the option
% `'deviation='` false (default). When reporting the results, plot both the
% simulated shock against the steady-state (balanced-growth path) database:
% The `&` operator combines two databases so that every time series has two
% columns.

d = sstatedb(m,startDate:endDate);
d.Ey(startDate) = log(1.01);
s = simulate(m,d,1:40);
s = dbextend(d,s);

qplot('simulate_simple_shock.q',d & s,plotrng,'tight',true);
grfun.ftitle('Consumption Demand Shock -- Full Levels');

%% Help on IRIS Functions Used in This File
%
% Use either `help` to display help in the command window, or `idoc`
% to display help in an HTML browser window.
%
%    help model/dbextend
%    help model/simulate
%    help model/sstatedb
%    help model/zerodb
%    help qreport/qplot
%    help grfun/ftitle
%    help dbplot
%    help dbextend
